In [1]:
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import quandl
import datetime
from scipy.stats import norm
import math

plt.rcParams['figure.figsize'] = (15,8)
%matplotlib inline
In [2]:
## This is small code to correct btc.csv file from coindesk. Because that csv has one extra comma at end.

#f = open('btc2.csv','a')
#for line in open('btc.csv').readlines():
#    l = line.strip('\n').strip(',')
#    f.write(l+'\n')
#f.close()
In [3]:
df = pd.read_csv('bt.csv',index_col='Date', parse_dates=True,
                 engine='python')
df2 = pd.read_csv('BTC.csv', index_col='date',parse_dates=True,names=['date','ranknow','open','high','low','close'],
                 engine='python',skiprows=1)
df3 = pd.read_csv('btc2.csv',index_col='date',parse_dates=True,
                  names=['date','txVolume(USD)','adjustedTxVolume(USD)','txCount','marketcap(USD)','price(USD)',
                         'exchangeVolume(USD)', 'generatedCoins','fees','activeAddresses', 'averageDifficulty',
                         'paymentCount','medianTxValue(USD)', 'medianFee','blockSize','blockCount'],
                  skiprows=1)
In [4]:
df3.dtypes
Out[4]:
txVolume(USD)            float64
adjustedTxVolume(USD)    float64
txCount                    int64
marketcap(USD)           float64
price(USD)               float64
exchangeVolume(USD)      float64
generatedCoins           float64
fees                     float64
activeAddresses            int64
averageDifficulty        float64
paymentCount             float64
medianTxValue(USD)       float64
medianFee                float64
blockSize                  int64
blockCount                 int64
dtype: object
In [5]:
df3.dropna(subset=['txVolume(USD)','marketcap(USD)'],inplace=True)
df3.fillna(0,inplace=True)
In [6]:
df3.dtypes
Out[6]:
txVolume(USD)            float64
adjustedTxVolume(USD)    float64
txCount                    int64
marketcap(USD)           float64
price(USD)               float64
exchangeVolume(USD)      float64
generatedCoins           float64
fees                     float64
activeAddresses            int64
averageDifficulty        float64
paymentCount             float64
medianTxValue(USD)       float64
medianFee                float64
blockSize                  int64
blockCount                 int64
dtype: object
In [13]:
df3[['price(USD)','txVolume(USD)','marketcap(USD)']].describe()
Out[13]:
price(USD) txVolume(USD) marketcap(USD)
count 1928.000000 1.928000e+03 1.928000e+03
mean 2048.630711 2.671062e+09 3.362301e+10
std 3389.588819 5.287822e+09 5.749219e+10
min 68.500000 3.287234e+07 7.792550e+08
25% 291.902500 2.430015e+08 4.151355e+09
50% 577.670000 7.047429e+08 7.671265e+09
75% 1222.840000 2.656065e+09 1.982230e+10
max 19475.800000 4.835307e+10 3.261413e+11

Bitcoin Closing prices plot (X = timespan, Y = log(ClosingPrice))

In [8]:
df.plot(figsize=(14,7),logy=True)
Out[8]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd3271aa9b0>

Bitcoin Open, High, Low and Close Price Plot.

In [9]:
df2.plot(figsize=(14,8),logy=True)
/home/sunny/anaconda2/envs/p3/lib/python3.6/site-packages/matplotlib/axes/_base.py:3443: UserWarning: Attempting to set identical bottom==top results
in singular transformations; automatically expanding.
bottom=1.0, top=1.0
  'bottom=%s, top=%s') % (bottom, top))
Out[9]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2e2697e10>
In [10]:
#bitcoin = quandl.get("BCHAIN/MKPRU")
In [11]:
#bitcoin['Value'].replace(0, np.nan, inplace=True)
#bitcoin = bitcoin.dropna()
#bitcoin.Value.plot(logy=True,figsize=(15,8));
#plt.legend()

Calculation of Returns from Closing Price for monte carlo simulation

In [12]:
df['Returns'] = (df['Daily Closing Price'].pct_change()+1).fillna(0)
#df.head()
dates = pd.date_range(start=datetime.date.today(),end='2018-12-30')
np.random.seed(1234)
simulated_returns_bitcoin = np.random.choice(df.Returns, size=(len(dates), 100))
sim_bitcoin_returns = pd.DataFrame(data=simulated_returns_bitcoin, index=dates)
#print(sim_bitcoin_returns.shape)
cum_sim_bitcoin = sim_bitcoin_returns.cumprod(axis=0)
#print(cum_sim_bitcoin.shape)
future = pd.DataFrame(data=cum_sim_bitcoin, index=dates)
future = future * df['Daily Closing Price'][-1]
future.head()
Out[12]:
0 1 2 3 4 5 6 7 8 9 ... 90 91 92 93 94 95 96 97 98 99
2018-08-12 7586.932163 7519.682051 7189.041223 7316.908367 7438.958927 7037.286226 8839.626241 8162.073011 7331.690000 7390.689338 ... 7284.109711 7448.894069 7323.448627 7452.853276 7606.526527 7140.880113 7159.696958 7667.065641 7506.617672 7678.185337
2018-08-13 6929.963177 7557.334250 7818.060955 7286.673208 7442.575338 7262.797755 8697.100639 8072.625636 7218.778171 7287.283006 ... 7480.701681 6307.655640 7362.298109 7541.842569 8035.364769 6793.083949 7293.329118 7456.387448 8171.761010 7404.916945
2018-08-14 6990.424404 7585.444879 8127.229296 7531.163142 7806.257462 7218.064739 8684.235106 8923.643479 7298.947143 7115.817523 ... 7604.463618 6294.206481 7135.765859 7414.915277 7992.923238 6729.434951 7371.101377 7747.217965 8272.045495 7404.916945
2018-08-15 7021.378313 7620.397011 8283.842812 7500.423701 7723.564057 7365.645989 8212.586826 9704.435750 7669.572951 7340.988149 ... 7462.324112 6215.835038 6760.489696 7589.857134 8407.961237 6664.632985 6557.114722 7673.434937 8285.538545 7404.916945
2018-08-16 7092.301326 7107.121024 8446.271103 7463.887190 6781.666001 7448.909667 8212.586826 9743.507022 7585.013161 6973.938742 ... 7448.619200 6227.620019 6807.524049 7495.489461 8567.001573 6720.507877 5838.661014 7547.923103 8270.391309 7480.576559

5 rows × 100 columns

Monte Carlo Price simulation plot for Returns

In [13]:
cum_sim_bitcoin.plot(logy=True,legend=False,figsize=(15,8))
Out[13]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2e207ff98>

Monte Carlo Simulation for Closing Price

In [14]:
future.plot(legend=False, figsize=(15,8))
Out[14]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2e22d7358>

Monte Carlo Simulation with log of price

In [15]:
future.plot(logy=True,legend=False, figsize=(15,8),
            title='Bitcoin price Monte Carlo simulations until Dec 31st 2018')
plt.xlabel('Date')
plt.ylabel("Price ($)")
Out[15]:
Text(0,0.5,'Price ($)')

Timeseries Graph of Transaction Volume and Network Value (Market Cap)

In [16]:
df3[['txVolume(USD)','marketcap(USD)']].plot(figsize=(15,8),logy=True)
Out[16]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2e1d435f8>

Timeseries graph of Exchange Volume and Market Cap.

In [17]:
df3[['exchangeVolume(USD)','marketcap(USD)']].plot(figsize=(15,8),logy=True)
Out[17]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2e0050390>

Calculation of NVT Ratio

In [18]:
df3['NVT'] = df3['marketcap(USD)'] / df3['txVolume(USD)']
df3.fillna(0,inplace=True)
df3.tail(5)['NVT']
Out[18]:
date
2018-08-03    28.554939
2018-08-04    35.576221
2018-08-05    35.103984
2018-08-06    28.528188
2018-08-07    26.307649
Name: NVT, dtype: float64

Timeseries plot of NVT.

In [19]:
df3[['NVT']].plot(figsize=(15,8))
Out[19]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2dbeb6e80>

Smothing NVT with 14 days rolling average

In [20]:
o = np.convolve(list(reversed(list(df3.NVT.values))), np.ones((14,))/14., mode='same')
o = list(reversed(o))
#print(o)
df3['SMOOTH_NVT'] = o
df3['SMOOTH_NVT'] = np.convolve(df3.NVT.values, np.ones((14,))/14., mode='same')
df3['SMOOTH_NVT'][:5]
Out[20]:
date
2013-04-28     8.906541
2013-04-29    10.476616
2013-04-30    11.551713
2013-05-01    12.473677
2013-05-02    13.563958
Name: SMOOTH_NVT, dtype: float64

Timeseries plot of NVT and Smooth NVT

In [21]:
plt.figure(figsize=(15,8))
plt.plot(df3.index,df3.NVT,color='grey',label='NVT')
plt.plot(df3.index, df3.SMOOTH_NVT,color = 'red',label='SMOOTH NVT')
plt.legend(loc='best')
plt.xlabel('Date')
plt.ylabel('NVT/SMOOTH NVT Ratio')
plt.title('NVT/SMOOTH NVT Ratio')
#df3[['NVT','SMOOTH_NVT']].plot(figsize=(15,8),colormap=plt.cm.Set2_r)
Out[21]:
Text(0.5,1,'NVT/SMOOTH NVT Ratio')

Time Series Plots of NVT, SMOOTH NVT and Scaled Market Cap & Log version as well.

In [22]:
df3['MCAP'] = df3['marketcap(USD)'] / 10e9
plt.figure(figsize=(15,8))
plt.plot(df3.index,df3.NVT,color='grey',label='NVT')
plt.plot(df3.index,df3.SMOOTH_NVT,color='black',label='SMOOTH NVT')
plt.plot(df3.index,df3.MCAP,color='orange',label='Scaled Market Cap(in Bn)')
plt.legend(loc='best')
plt.xlabel('Date')
plt.ylabel('NVT/SMOOTH NVT/ Market Cap')
df3[['NVT','SMOOTH_NVT','MCAP']].plot(figsize=(15,8),colormap=plt.cm.Set2_r,logy=True)
Out[22]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2e00f0780>

Timeseries plot of Exchange Volume (USD) and Market Cap (USD)

In [23]:
df3.index.name = 'Date'
df3[['exchangeVolume(USD)','marketcap(USD)']].plot(figsize=(15,8))
Out[23]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2dbc33128>

Timeseries plot of Scaled Exchange Volume (USD) and Market Cap (USD)

In [24]:
plt.figure(figsize=(15,8))
plt.plot(df3.index,df3['exchangeVolume(USD)']/1e9,color='orange',label='Exchange Volume (in Bn)')
plt.plot(df3.index,df3['marketcap(USD)']/1e9,color='blue',label='Market Cap (in Bn)')
plt.legend(loc='best')
plt.xlabel('Date')
plt.ylabel('Exchange Vol/ Market Cap')
df3[['exchangeVolume(USD)','marketcap(USD)']].plot(figsize=(15,8))
Out[24]:
<matplotlib.axes._subplots.AxesSubplot at 0x7fd2dbad5be0>
In [25]:
#df3.reset_index(inplace=True)
#df3.Date = [datetime.date(d.year,d.month,1) for d in df3.Date]
#df3 = df3.groupby(['Date']).sum()
#df3[['NVT','SMTH_NVT','MCAP']].plot(figsize=(15,8),colormap=plt.cm.Set2_r)

Timeseries plot of Transaction Volume (USD)

In [26]:
plt.figure(figsize=(12,6))
plt.plot(df3.index.values,df3['txVolume(USD)'].values)
l,t = plt.yticks([0,10e9,20e9,30e9,40e9,50e9],['0 Bn','10 Bn','20 Bn','30 Bn','40 Bn','50 Bn'])                 
#plt.plot(df3.index.values,df3['marketcap(USD)'].values)
plt.xlabel('Year')
plt.ylabel('txVolume(USD)')
#plt.axes()
#plt.plot(df3.index.values,df3['marketcap(USD)'].values)
Out[26]:
Text(0,0.5,'txVolume(USD)')

Simulating NVT Ratio using Monte Carlo Simulation

In [27]:
df3['NVT_R'] = (df3['NVT'].pct_change()+1).fillna(0)
dates = pd.date_range(start=datetime.date.today(),end='2018-12-30')
np.random.seed(1234)
simulated_nvt = np.random.choice(df3.NVT_R, size=(len(dates), 100))
sim_nvt = pd.DataFrame(data=simulated_nvt, index=dates)
cum_sim_nvt = sim_nvt.cumprod(axis=0)
future2 = pd.DataFrame(data=cum_sim_nvt, index=dates)
future2 = future2 * df3['NVT'][-1]
In [28]:
future2.plot(legend=False, figsize=(15,8),
            title='Bitcoin NVT ratio Monte Carlo simulations until Dec 31st 2018')
plt.xlabel('Date')
plt.ylabel("NVT")
Out[28]:
Text(0,0.5,'NVT')

Simulating Transaction Volume and Network Value(Market Cap) using Monte Carlo

In [29]:
df3['txVolume_R'] = (df3['txVolume(USD)'].pct_change()+1).fillna(0)
dates = pd.date_range(start=datetime.date.today(),end='2018-12-30')
np.random.seed(1234)
simulated_tx_vol = np.random.choice(df3.txVolume_R, size=(len(dates), 100))
sim_tx_vol = pd.DataFrame(data=simulated_tx_vol, index=dates)
cum_sim_tx_vol = sim_tx_vol.cumprod(axis=0)
future3 = pd.DataFrame(data=cum_sim_tx_vol, index=dates)
future3 = future3 * df3['txVolume(USD)'][-1]

df3['marketcap_R'] = (df3['marketcap(USD)'].pct_change()+1).fillna(0)
dates = pd.date_range(start=datetime.date.today(),end='2018-12-30')
np.random.seed(1234)
simulated_marketcap_R = np.random.choice(df3.marketcap_R, size=(len(dates), 100))
sim_marketcap_R = pd.DataFrame(data=simulated_marketcap_R, index=dates)
cum_sim_marketcap_R = sim_marketcap_R.cumprod(axis=0)
future4 = pd.DataFrame(data=cum_sim_marketcap_R, index=dates)
future4 = future4 * df3['marketcap(USD)'][-1]

df4 = pd.DataFrame(data = future4.values / future3.values)
df4 = df4.set_index(future4.index)
/home/sunny/anaconda2/envs/p3/lib/python3.6/site-packages/ipykernel_launcher.py:19: RuntimeWarning: invalid value encountered in true_divide
In [30]:
future3.plot(legend=False, figsize=(15,8), title='Bitcoin Transaction Volume Monte Carlo simulations until Dec 31st 2018')
plt.xlabel('Date')
plt.ylabel("Transaction Volume")
future4.plot(legend=False, figsize=(15,8), title='Bitcoin market Cap Monte Carlo simulations until Dec 31st 2018')
plt.xlabel('Date')
plt.ylabel("Market Cap")
df4.plot(legend=False, figsize=(15,8),
            title='Bitcoin NVT ratio Monte Carlo simulations until Dec 31st 2018')
plt.Circle((2000,'2018-09-01'),20,visible=True,fill=True)
plt.xlabel('Date')
plt.ylabel("NVT")
Out[30]:
Text(0,0.5,'NVT')
In [31]:
df3.head()
Out[31]:
txVolume(USD) adjustedTxVolume(USD) txCount marketcap(USD) price(USD) exchangeVolume(USD) generatedCoins fees activeAddresses averageDifficulty ... medianTxValue(USD) medianFee blockSize blockCount NVT SMOOTH_NVT MCAP NVT_R txVolume_R marketcap_R
Date
2013-04-28 6.879868e+07 3.153714e+07 41702 1.500520e+09 135.30 0.0 4425.0000 32.791942 117984 8.974296e+06 ... 49.980091 0.0005 21597536 177 21.810300 8.906541 0.150052 0.000000 0.000000 0.000000
2013-04-29 1.138128e+08 4.821652e+07 51602 1.491160e+09 134.44 0.0 4349.9988 45.724114 86925 8.974296e+06 ... 62.425891 0.0005 25676453 174 13.101861 10.476616 0.149116 0.600719 1.654288 0.993762
2013-04-30 8.426632e+07 5.686790e+07 47450 1.597780e+09 144.00 0.0 3725.0000 45.748651 76871 9.854414e+06 ... 26.942145 0.0005 27042465 149 18.961074 11.551713 0.159778 1.447205 0.740394 1.071501
2013-05-01 1.206825e+08 8.208736e+07 55176 1.542820e+09 139.00 0.0 3775.0000 40.885938 83564 1.007629e+07 ... 38.671980 0.0005 25727393 151 12.784120 12.473677 0.154282 0.674230 1.432156 0.965602
2013-05-02 9.337533e+07 5.823736e+07 55295 1.292190e+09 116.38 0.0 3350.0000 52.554004 81920 1.007629e+07 ... 35.382348 0.0005 26388234 134 13.838666 13.563958 0.129219 1.082489 0.773727 0.837551

5 rows × 21 columns

In [32]:
future2.plot(legend=False, figsize=(15,8),
            title='Bitcoin NVT ratio Monte Carlo simulations until Dec 31st 2018')
plt.xlabel('Date')
plt.ylabel("NVT")
#plt.axhline(y=50000,xmin=0,xmax=10000,color='blue')
#plt.axvline(x=50,ymin=0,ymax=10000,color='blue')
Out[32]:
Text(0,0.5,'NVT')

Functoin plots of e^−2αXt with α = 0.1, 0.5, 1, 2. Green zone: positive drift; red zone:negative drift.

In [33]:
fig = plt.figure(figsize=(8,6))
Xt = np.linspace(0,2.5,100)
c = 0.5
for alpha in [0.1, 0.5, 1.0, 2.0]:
    plt.plot(Xt,np.exp(Xt*alpha*-2),label='Alpha %s'%alpha)    
t = zip(Xt,np.exp(Xt*alpha*-2))
location = [(j,i) for j,i in enumerate(t) if  str(i[1])[:3] == '0.5']
#print(location)
plt.plot(Xt[location[0][0]:],[c]*len(Xt[location[0][0]:]),color='black',label='c %s'%c)

y2 = np.array(len(Xt)*[c])
upper = zip(Xt,np.exp(Xt*2.1*-2))

plt.fill_between(Xt, y2 ,np.exp(Xt*0.1*-2), alpha=0.3)
plt.fill_between(Xt, y2, np.exp(Xt*2.0*-2), where= Xt > location[0][1][0], alpha=0.3)
plt.fill_between(Xt, np.exp(Xt*2.0*-2), 0 , color='white')
plt.xlabel('Xt')
plt.ylabel('e^-2aXt')
plt.legend(loc='best')
Out[33]:
<matplotlib.legend.Legend at 0x7fd2db5172b0>

Sample path for Xt in 4 years time. Parameters are chosen as α = 1, X0 = 0 and dt 1/250 .

In [57]:
def get_brownian(n):
    # Process parameters
    delta = 0.25
    #delta = 1.0/250.0
    dt = 1.0/250.0
    #dt = 1.0
    # Initial condition.
    x = 0.0
    Wt = [0.0]
    # Iterate to compute the steps of the Brownian motion.
    for k in range(n):
        x = x + norm.rvs(scale=delta**2*dt)
        #x = norm.rvs(scale=delta**2*dt)
        Wt.append(x)
    return Wt

def brownian_path(N):
    Δt_sqrt = math.sqrt(1 / N)
    Z = np.random.randn(N)
    Z[0] = 0
    B = np.cumsum(Δt_sqrt * Z)
    return B

def brownian(x0, n, dt, delta, out=None):
    """
    Generate an instance of Brownian motion (i.e. the Wiener process):

        X(t) = X(0) + N(0, delta**2 * t; 0, t)

    where N(a,b; t0, t1) is a normally distributed random variable with mean a and
    variance b.  The parameters t0 and t1 make explicit the statistical
    independence of N on different time intervals; that is, if [t0, t1) and
    [t2, t3) are disjoint intervals, then N(a, b; t0, t1) and N(a, b; t2, t3)
    are independent.
    
    Written as an iteration scheme,

        X(t + dt) = X(t) + N(0, delta**2 * dt; t, t+dt)


    If `x0` is an array (or array-like), each value in `x0` is treated as
    an initial condition, and the value returned is a numpy array with one
    more dimension than `x0`.

    Arguments
    ---------
    x0 : float or numpy array (or something that can be converted to a numpy array
         using numpy.asarray(x0)).
        The initial condition(s) (i.e. position(s)) of the Brownian motion.
    n : int
        The number of steps to take.
    dt : float
        The time step.
    delta : float
        delta determines the "speed" of the Brownian motion.  The random variable
        of the position at time t, X(t), has a normal distribution whose mean is
        the position at time t=0 and whose variance is delta**2*t.
    out : numpy array or None
        If `out` is not None, it specifies the array in which to put the
        result.  If `out` is None, a new numpy array is created and returned.

    Returns
    -------
    A numpy array of floats with shape `x0.shape + (n,)`.
    
    Note that the initial value `x0` is not included in the returned array.
    """

    x0 = np.asarray(x0)

    # For each element of x0, generate a sample of n numbers from a
    # normal distribution.
    print(x0.shape + (n,))
    r = norm.rvs(size=x0.shape + (n,), scale=delta*math.sqrt(dt))

    # If `out` was not given, create an output array.
    if out is None:
        out = np.empty(r.shape)

    # This computes the Brownian motion by forming the cumulative sum of
    # the random samples. 
    np.cumsum(r, axis=-1, out=out)

    # Add the initial condition.
    out += np.expand_dims(x0, axis=-1)

    return out
In [66]:
def get_series():
    dt = 1.0 / 250.0
    tim = np.linspace(0.0,4.0,4000)
    alpha = 1.0
    eta = 0.1
    c = 0.5
    Xt = np.zeros(len(tim)+1)
    dWt = brownian(0, len(tim), dt, 0.5)
    #print(dWt)
    eq,eq10 = np.log(0.5) / (-2.0 * alpha),np.log(0.05) / (-2.0 * alpha)
    for i in range(len(tim)):
        t = np.exp(-2.0 * alpha * Xt[i])
        tf = t - c
        dXt = (eta * tf * dt) + dWt[i]
        Xt[i+1] = dXt
    return np.array(Xt),tim,eq,eq10
In [69]:
#Xt,eq,eq10 = get_series()
Xt,tim,eq,eq10 = get_series()
#Xt = np.exp(Xt)
fig = plt.figure(figsize=(8,6))
plt.plot(tim,Xt[:len(tim)])
plt.plot(tim, len(tim)*[0.0],linestyle='dashed',label='Intial Point',color='green')

plt.plot(tim, len(tim)*[eq10],linestyle='dashed',label = '10% of Equilbrium Level',color='red')

plt.plot(tim, len(tim)*[eq],linestyle='dashed',label = 'Equilbrium Level',color='black')
plt.xlabel('Time (Years)')
plt.ylabel('Xt')
plt.legend(loc='best')
(4000,)
Out[69]:
<matplotlib.legend.Legend at 0x7fd2da8e3eb8>
In [72]:
ixic = pd.read_csv('^IXIC.csv',index_col='Date',parse_dates=True)
ssec = pd.read_csv('^SSEC.csv',index_col='Date',parse_dates=True)
In [73]:
ixic.head(5)
Out[73]:
Open High Low Close Adj Close Volume
Date
1971-02-05 100.000000 100.000000 100.000000 100.000000 100.000000 0
1971-02-08 100.839996 100.839996 100.839996 100.839996 100.839996 0
1971-02-09 100.760002 100.760002 100.760002 100.760002 100.760002 0
1971-02-10 100.690002 100.690002 100.690002 100.690002 100.690002 0
1971-02-11 101.449997 101.449997 101.449997 101.449997 101.449997 0
In [74]:
ssec.head(5)
Out[74]:
Open High Low Close Adj Close Volume
Date
1997-07-02 1255.909058 1261.571045 1147.331055 1199.061035 1199.061035 0.0
1997-07-03 1194.676025 1194.676025 1149.939941 1150.623047 1150.623047 0.0
1997-07-04 1138.921021 1163.249023 1124.776001 1159.342041 1159.342041 0.0
1997-07-07 1161.707031 1163.447021 1085.572021 1096.818970 1096.818970 0.0
1997-07-08 1092.798950 1115.432983 1066.043945 1109.666016 1109.666016 0.0

NASDAQ ^IXIC Calculation and Graphs

In [253]:
t1 = [int(i) for i in '1997-01-02'.split('-')]   # (P0 = 1, 280), (Pt = 1, 436)
t2 = [int(i) for i in '1997-06-26'.split('-')]

X1_i = ixic[(ixic.index > pd.to_datetime('1997-01-02')) & (ixic.index < pd.to_datetime('1997-06-26'))] 

t3 = [int(i) for i in '1997-06-26'.split('-')] # (Pt = 1, 436), (P t = 2, 309)
t4 = [int(i) for i in '1999-02-10'.split('-')]

X2_i = ixic[(ixic.index > pd.to_datetime('1997-06-26')) & (ixic.index < pd.to_datetime('1999-02-10'))] 

t5 = [int(i) for i in '1999-02-10'.split('-')]  #(P = 2, 309), (P = 3, 171).
t6 = [int(i) for i in '2000-10-18'.split('-')]

X3_i = ixic[(ixic.index > pd.to_datetime('1999-02-10')) & (ixic.index < pd.to_datetime('2000-10-18'))] 

XR_hat_i = 0.67 # (PR = 2, 502).
eta_hat_i = 0.39
alpha_hat_i = 0.23
sigma_hat_i = 0.43
c_hat_i = 0.73

##  OU Process
k_hat_i, mu_hat_i, sig_hat_i = (0.47, 1.09, 0.31)

## BM Process
mw_hat_i, sigm_hat_i = (0.25, 0.31)

total_i = (len(X1_i)+len(X2_i)+len(X3_i))
dt_i = 1.0 / total_i

adj_price_i = np.log(np.array(list(X1_i['Adj Close']/1000)+list(X2_i['Adj Close']/1000)+list(X3_i['Adj Close']/1000)))
timespan_i = list(X1_i.index)+list(X2_i.index)+list(X3_i.index)
In [303]:
def get_Xt_path(X0, dt , eta_hat,alpha_hat,sigma_hat,c_hat, total ):
    tim = np.arange(0.0,total,1)
    Xt = np.zeros(len(tim))
    dWt = brownian_path(len(tim))
    Xt[0] = X0
    for i in range(len(tim)-1):
        t = np.exp(-2.0 * alpha_hat * (Xt[i]/sigma_hat)) - c_hat
        dXt = (eta_hat * t * dt) + dWt[i]
        Xt[i+1] = dXt
    return np.array(Xt),tim

def get_Bm_path(X0, dt , mw_hat, sigma_hat, total):
    #t = np.linspace(0, T, N)
    W = np.random.standard_normal(size = total) 
    W = np.cumsum(W)*np.sqrt(dt) ### standard brownian motion ###
    X = (mw_hat-0.5*sigma_hat**2)*dt + sigma_hat*W 
    X = X0*np.exp(X)
    return X

def get_Ou_path(dt, theta, mu_hat,sig_hat, total ):
    tim = np.arange(0.0,total,1)
    Xt = np.zeros(len(tim))
    dWt = brownian_path(len(tim))
    for i in range(len(tim)-1):
        dXt = theta * (mu_hat - Xt[i]) *dt + sig_hat * dWt[i]
        Xt[i+1] = dXt
    return np.array(Xt)
In [317]:
Ou_i = get_Ou_path(dt_i, k_hat_i, mu_hat_i ,sig_hat_i,total_i)
plt.figure(figsize = (8,6))
plt.plot(timespan_i,Ou_i,color='orange',label='Ornstein Ulhenback Samplt Path',linestyle='--')
plt.legend(loc='best')
Out[317]:
<matplotlib.legend.Legend at 0x7fd2d2694710>
In [276]:
Xt_i,tim_i = get_Xt_path(0, dt_i,eta_hat_i,alpha_hat_i,sigma_hat_i,c_hat_i,total_i)
plt.figure(figsize = (8,6))
plt.plot(timespan_i,Xt_i,color='blue',label='Xt Sample path',linestyle='--')
plt.legend(loc='best')
Out[276]:
<matplotlib.legend.Legend at 0x7fd2d33b72b0>
In [294]:
bm_i = get_Bm_path(XR_hat_i,dt,mw_hat_i,sigm_hat_i, total_i)
plt.figure(figsize = (8,6))
plt.plot(timespan_i,np.log(bm_i),color='green',label='Brownian Motion Sample path',linestyle='--')
plt.legend(loc='best')
Out[294]:
<matplotlib.legend.Legend at 0x7fd2d2cb1e10>
In [318]:
plt.figure(figsize = (8,6))

plt.plot(timespan_i,Xt_i,color='blue',label='Xt Sample path',linestyle='--')
plt.plot(timespan_i,np.log(bm_i),color='green',label='Brownian Motion Sample path',linestyle='--')
plt.plot(timespan_i, adj_price_i,color='red',label='Log of adj close price',linestyle='--')
plt.plot(timespan_i, Ou_i,color='orange',label='Ornstein-Uhlenback Sample path',linestyle='--')
plt.legend(loc='best')
Out[318]:
<matplotlib.legend.Legend at 0x7fd2d261fa58>
In [238]:
plt.figure(figsize=(15,8))
plt.plot(X1.index, np.log(X1_i.Close),X2_i.index, np.log(X2_i.Close),X3_i.index, np.log(X3_i.Close),color='blue')
plt.fill_between(X1_i.index,np.max(np.log(X1_i.Close)),np.min(np.log(X1_i.Close)),alpha=0.4,color='green')
plt.fill_between(X2_i.index,np.max(np.log(X2_i.Close)),np.min(np.log(X2_i.Close)),alpha=0.6,color='yellow')
plt.fill_between(X3_i.index,np.max(np.log(X3_i.Close)),np.min(np.log(X3_i.Close)),alpha=0.6,color='red')

y = pd.date_range(start='2000-10-18',end='2001-10-10')
X4_i = ixic[ixic.index.isin(y)]
for i in range(1000):
    Xt,tim = get_series2(0,dt,eta_hat,alpha_hat,sigma_hat,c_hat,total)
    plt.plot(X4_i.index, Xt[:len(X4_i.index)]+np.log(X3_i.Close[-1]),color='grey')
plt.plot(X4_i.index, np.log(X4_i.Close),color='red')
Out[238]:
[<matplotlib.lines.Line2D at 0x7fd2d430d358>]

Shangai calculation and Graphs

In [320]:
# Shangai 
t1 = [int(i) for i in '2006-01-04'.split('-')]   # (P0 = 1, 180), (P t = 1,288)
t2 = [int(i) for i in '2006-03-06'.split('-')]

X1_s = ssec[(ssec.index > pd.to_datetime('2006-01-04')) & (ssec.index < pd.to_datetime('2006-03-06'))] 

t3 = [int(i) for i in '2006-03-06'.split('-')] # (Pt1 = 1, 288), (Pt2 = 4, 053)
t4 = [int(i) for i in '2007-05-30'.split('-')]

X2_s = ssec[(ssec.index > pd.to_datetime('2006-03-06')) & (ssec.index < pd.to_datetime('2007-05-30'))] 

t5 = [int(i) for i in '2007-05-30'.split('-')]  # (P = 4, 053) to  (P = 3, 116).
t6 = [int(i) for i in '2008-04-21'.split('-')]

X3_s = ssec[(ssec.index > pd.to_datetime('2007-05-30')) & (ssec.index < pd.to_datetime('2008-04-21'))] 

XR_hat_s = 1.23 # (P R = 4, 040).
eta_hat_s = 0.32
alpha_hat_s = 0.14
sigma_hat_s = 0.56
c_hat_s = 0.70

##  OU Process
k_hat_s, mu_hat_s,sig_hat_s = (3.30, 0.97, 1.20)

## BM Process
mw_hat_s, sigm_hat_s = (0.44, 0.33)

total_s = (len(X1_s)+len(X2_s)+len(X3_s))
dt_s = 1.0 / total_s
timespan_s = list(X1_s.index)+list(X2_s.index)+list(X3_s.index)
adj_price_s = np.log(np.array(list(X1_s['Adj Close']/1000)+list(X2_s['Adj Close']/1000)+list(X3_s['Adj Close']/1000)))
In [321]:
Ou_s = get_Ou_path(dt_s, k_hat_s, mu_hat_s ,sig_hat_s,total_s)
plt.figure(figsize = (8,6))
plt.plot(timespan_s, Ou_s, color='orange',label='Ornstein Ulhenback Samplt Path',linestyle='--')
plt.legend(loc='best')
Out[321]:
<matplotlib.legend.Legend at 0x7fd2d2587470>
In [334]:
Xt_s,tim_s = get_Xt_path(0, dt_s,eta_hat_s,alpha_hat_s,sigma_hat_s,c_hat_s,total_s)
plt.figure(figsize = (8,6))
plt.plot(timespan_s,Xt_s,color='blue',label='Xt Sample path',linestyle='--')
plt.legend(loc='best')
Out[334]:
<matplotlib.legend.Legend at 0x7fd2d204e2b0>
In [338]:
bm_s = get_Bm_path(XR_hat_s,dt,mw_hat_s,sigm_hat_s, total_s)
plt.figure(figsize = (8,6))
plt.plot(timespan_s, np.log(bm_s),color='green',label='Brownian Motion Sample path',linestyle='--')
plt.legend(loc='best')
Out[338]:
<matplotlib.legend.Legend at 0x7fd2d1f795f8>
In [339]:
plt.figure(figsize = (8,6))

plt.plot(timespan_s,Xt_s,color='blue',label='Xt Sample path',linestyle='--')
plt.plot(timespan_s,np.log(bm_s),color='green',label='Brownian Motion Sample path',linestyle='--')
plt.plot(timespan_s, adj_price,color='red',label='Log of adj close price',linestyle='--')
plt.plot(timespan_s, Ou_s,color='orange',label='Ornstein-Uhlenback Sample path',linestyle='--')
plt.legend(loc='best')
Out[339]:
<matplotlib.legend.Legend at 0x7fd2d1dfdf98>
In [212]:
plt.figure(figsize=(15,8))
plt.plot(X1.index, np.log(X1.Close),X2.index, np.log(X2.Close),X3.index, np.log(X3.Close),color='blue')
plt.fill_between(X1.index,np.max(np.log(X1.Close)),np.min(np.log(X1.Close)),alpha=0.4,color='green')
plt.fill_between(X2.index,np.max(np.log(X2.Close)),np.min(np.log(X2.Close)),alpha=0.6,color='yellow')
plt.fill_between(X3.index,np.max(np.log(X3.Close)),np.min(np.log(X3.Close)),alpha=0.6,color='red')

y = pd.date_range(start='2008-04-18',end='2008-11-21')
X4 = ssec[ssec.index.isin(y)]
for i in range(1000):
    Xt,tim = get_series2(0,dt,eta_hat,alpha_hat,sigma_hat,c_hat,total)
    plt.plot(X4.index, Xt[:len(X4.index)]+np.log(X3.Close[-1]),color='grey')
plt.plot(X4.index, np.log(X4.Close),color='red')
Out[212]:
[<matplotlib.lines.Line2D at 0x7f6937ba6e48>]

Bitcoin Xt path calculation and Graphs

In [223]:
t1 = [int(i) for i in '2016-01-01'.split('-')]
t2 = [int(i) for i in '2016-05-30'.split('-')]

X1_b = df[(df.index.date > datetime.date(*t1)) & (df.index.date < datetime.date(*t2))] #(P0 = 433), (Pt1 = 528)

t3 = [int(i) for i in '2016-05-30'.split('-')]
t4 = [int(i) for i in '2017-08-13'.split('-')]

X2_b = df[(df.index.date > datetime.date(*t3)) & (df.index.date < datetime.date(*t4))] #(Pt = 528), (Pt = 4,327);

t5 = [int(i) for i in '2017-08-13'.split('-')]
t6 = [int(i) for i in '2017-12-10'.split('-')]

X3_b = df[(df.index.date > datetime.date(*t5)) & (df.index.date < datetime.date(*t6))] # (P = 4,327), (P = 14,371).

XR_hat_b = 2.30 # (PR = 4, 327);
eta_hat_b = 0.51
alpha_hat_b = 0.08
sigma_hat_b = 0.91
c_hat_b = 0.69

total_b = (len(X1_b)+len(X2_b)+len(X3_b))
dt_b = 1.0 / total_b
In [227]:
Xt_b,tim_b = get_Xt_path(0, dt_b,eta_hat_b,alpha_hat_b,sigma_hat_b,c_hat_b,total_b)
#bm_b = get_Bm_path(XR_hat_b,dt,mu_hat_b,sigm_hat_b, total_b)

adj_price = np.log(np.array(list(X1_b['Daily Closing Price']/1000)+list(X2_b['Daily Closing Price']/1000)+list(X3_b['Daily Closing Price']/1000)))

plt.figure(figsize = (10,6))
timespan = list(X1_b.index)+list(X2_b.index)+list(X3_b.index)
plt.plot(timespan,Xt_b,color='blue',label='Xt Sample path',linestyle='--')
#plt.plot(timespan,bm_b,color='green',label='Brownian Motion Sample path',linestyle='--')
plt.plot(timespan, adj_price,color='red',label='Log of adj close price',linestyle='--')
#plt.plot(timespan, adj_price,color='orange',label='Ornstein-Uhlenback Sample path',linestyle='--')
plt.legend(loc='best')
Out[227]:
<matplotlib.legend.Legend at 0x7fd2d50182e8>
In [159]:
df['Daily Closing Price'].plot(figsize=(15,6))
plt.legend(loc='best')
Out[159]:
<matplotlib.legend.Legend at 0x7fd2d7be2d30>
In [231]:
#df3['txVolume(USD)'].plot(figsize=(15,6), kind='bar',color='red')
plt.figure(figsize=(15,6))
plt.bar(df3.index,df3['txVolume(USD)'].values,color='red',label='Volume')
plt.legend(loc='best')
Out[231]:
<matplotlib.legend.Legend at 0x7f6933d66a20>
In [361]:
plt.figure(figsize=(12,8))
dr = pd.date_range(start='2017-12-10',end='2018-01-12')
df_b = df[df.index.isin(dr)]
plt.plot(df_b.index,df_b['Daily Closing Price'],label='Historical Closing Price',color='red',linestyle='-.')
plt.plot(df_b.index,[df_b['Daily Closing Price'][0]]*len(df_b),label ='- 0.0%',linestyle='-',color='grey')
plt.plot(df_b.index,[13700]*len(df_b),label ='- 5.0%',linestyle='--',color='orange')
plt.plot(df_b.index,[12950]*len(df_b),label ='- 10.0%',linestyle='--',color='purple')
plt.plot(df_b.index,[12100]*len(df_b),label ='- 15.0%',linestyle='--',color='red')
plt.legend(loc='best')
Out[361]:
<matplotlib.legend.Legend at 0x7fd2d1448ac8>

Bitcoin Volatility Graph

In [349]:
bitcoin_volatility = pd.read_csv('bt.csv',index_col='Date', parse_dates=True,
                 engine='python')
bitcoin_volatility['Returns'] = bitcoin_volatility['Daily Closing Price'].pct_change().fillna(0)
bitcoin_volatility['Returns'] = bitcoin_volatility['Returns']*100

daily_volatility = bitcoin_volatility['Returns'].std()
sixty_day_volatility = daily_volatility * math.sqrt(6)
plt.figure(figsize=(15,8))
plt.plot(bitcoin_volatility.index,bitcoin_volatility.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,sixty_day_volatility))
plt.legend(loc='best')
Out[349]:
<matplotlib.legend.Legend at 0x7f2239d69f98>
In [348]:
bitcoin_volatility1 = bitcoin_volatility.reset_index()
bitcoin_volatility1.Date = [datetime.date(d.year,d.month,1) for d in bitcoin_volatility1.Date]
bitcoin_volatility1 = bitcoin_volatility1.set_index('Date')
bitcoin_volatility1 = bitcoin_volatility1[['Returns']]
bitcoin_volatility1 = bitcoin_volatility1.groupby(['Date']).std()
daily_volatility = bitcoin_volatility1.std()
sixty_day_volatility = daily_volatility * math.sqrt(6)

plt.figure(figsize=(15,8))
plt.plot(bitcoin_volatility1.index,bitcoin_volatility1.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,thirty_day_volatility))
plt.legend(loc='best')
Out[348]:
<matplotlib.legend.Legend at 0x7f2239e7c6a0>

Ripple Daily Volatility and 60 day volatility

In [350]:
ripple = pd.read_csv('ripple.csv',index_col='Date', parse_dates=True,
                 engine='python')
ripple['Returns'] = ripple['Close'].pct_change().fillna(0)
ripple['Returns'] = ripple['Returns']*100

daily_volatility = ripple['Returns'].std()
sixty_day_volatility = daily_volatility * math.sqrt(6)
plt.figure(figsize=(15,8))
plt.plot(ripple.index,ripple.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,sixty_day_volatility))
plt.legend(loc='best')
Out[350]:
<matplotlib.legend.Legend at 0x7f22400874e0>
In [352]:
ripple1 = ripple.reset_index()
ripple1.Date = [datetime.date(d.year,d.month,1) for d in ripple1.Date]
ripple1 = ripple1.set_index('Date')
ripple1 = ripple1[['Returns']]
ripple1 = ripple1.groupby(['Date']).std()
daily_volatility = ripple1.std()
sixty_day_volatility = daily_volatility * math.sqrt(12)

plt.figure(figsize=(15,8))
plt.plot(ripple1.index,ripple1.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,sixty_day_volatility))
plt.legend(loc='best')
Out[352]:
<matplotlib.legend.Legend at 0x7f2239d0afd0>

Lite coin Daily Volatility and 60 day volatility

In [353]:
litecoin = pd.read_csv('litecoin.csv',index_col='Date', parse_dates=True,
                 engine='python')
litecoin['Returns'] = litecoin['Close'].pct_change().fillna(0)
litecoin['Returns'] = litecoin['Returns']*100

daily_volatility = litecoin['Returns'].std()
sixty_day_volatility = daily_volatility * math.sqrt(6)
plt.figure(figsize=(15,8))
plt.plot(litecoin.index,litecoin.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,sixty_day_volatility))
plt.legend(loc='best')
Out[353]:
<matplotlib.legend.Legend at 0x7f224007f320>
In [354]:
litecoin1 = litecoin.reset_index()
litecoin1.Date = [datetime.date(d.year,d.month,1) for d in litecoin1.Date]
litecoin1 = litecoin1.set_index('Date')
litecoin1 = litecoin1[['Returns']]
litecoin1 = litecoin1.groupby(['Date']).std()
daily_volatility = litecoin1.std()
sixty_day_volatility = daily_volatility * math.sqrt(12)

plt.figure(figsize=(15,8))
plt.plot(litecoin1.index,litecoin1.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,sixty_day_volatility))
plt.legend(loc='best')
Out[354]:
<matplotlib.legend.Legend at 0x7f224040cac8>

Etherium Daily Volatility and 60 Day Volatility

In [355]:
etherium = pd.read_csv('etherium.csv',index_col='Date', parse_dates=True,
                 engine='python')
etherium['Returns'] = etherium['Close'].pct_change().fillna(0)
etherium['Returns'] = etherium['Returns']*100

daily_volatility = etherium['Returns'].std()
sixty_day_volatility = daily_volatility * math.sqrt(6)
plt.figure(figsize=(15,8))
plt.plot(etherium.index,etherium.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,sixty_day_volatility))
plt.legend(loc='best')
Out[355]:
<matplotlib.legend.Legend at 0x7f2239c2bc88>
In [356]:
etherium1 = etherium.reset_index()
etherium1.Date = [datetime.date(d.year,d.month,1) for d in etherium1.Date]
etherium1 = etherium1.set_index('Date')
etherium1 = etherium1[['Returns']]
etherium1 = etherium1.groupby(['Date']).std()
daily_volatility = etherium1.std()
sixty_day_volatility = daily_volatility * math.sqrt(12)

plt.figure(figsize=(15,8))
plt.plot(etherium1.index,etherium1.Returns,label='Daily Percentage change in price.')
plt.xlabel('Time')
plt.ylabel('Percentage change in price')
plt.title('Daily Volatility: %0.2f,  60 day volatility: %0.2f'%(daily_volatility,sixty_day_volatility))
plt.legend(loc='best')
Out[356]:
<matplotlib.legend.Legend at 0x7f2239c3f358>